{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Ground state solvers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "<img src=\"aux_files/H2_gs.png\" width=\"200\">\n",
    "\n",
    "In this tutorial we are going to discuss the ground state calculation interface of Qiskit Nature. The goal is to compute the ground state of a molecular Hamiltonian. This Hamiltonian can be electronic or vibrational. To know more about the preparation of the Hamiltonian, check out the Electronic structure and Vibrational structure tutorials. \n",
    "\n",
    "The first step is to define the molecular system. In the following we ask for the electronic part of a hydrogen molecule."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Couldn't find cython becke routine\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/manoel/opt/anaconda3/envs/Qiskitenv/lib/python3.8/site-packages/pyscf/lib/misc.py:46: H5pyDeprecationWarning: Using default_file_mode other than 'r' is deprecated. Pass the mode to h5py.File() instead.\n",
      "  h5py.get_config().default_file_mode = 'a'\n"
     ]
    }
   ],
   "source": [
    "from qiskit import Aer\n",
    "from qiskit_nature.drivers import UnitsType, Molecule\n",
    "from qiskit_nature.drivers.second_quantization import (\n",
    "    ElectronicStructureDriverType,\n",
    "    ElectronicStructureMoleculeDriver,\n",
    ")\n",
    "from qiskit_nature.problems.second_quantization import ElectronicStructureProblem\n",
    "from qiskit_nature.converters.second_quantization import QubitConverter\n",
    "from qiskit_nature.mappers.second_quantization import JordanWignerMapper\n",
    "\n",
    "molecule = Molecule(\n",
    "    geometry=[[\"H\", [0.0, 0.0, 0.0]], [\"H\", [0.0, 0.0, 0.735]]], charge=0, multiplicity=1\n",
    ")\n",
    "driver = ElectronicStructureMoleculeDriver(\n",
    "    molecule, basis=\"sto3g\", driver_type=ElectronicStructureDriverType.PYSCF\n",
    ")\n",
    "\n",
    "es_problem = ElectronicStructureProblem(driver)\n",
    "qubit_converter = QubitConverter(JordanWignerMapper())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The Solver\n",
    "\n",
    "Then we need to define a solver. The solver is the algorithm through which the ground state is computed. \n",
    "\n",
    "Let's first start with a purely classical example: the NumPy minimum eigensolver. This algorithm exactly diagonalizes the Hamiltonian. Although it scales badly, it can be used on small systems to check the results of the quantum algorithms. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from qiskit.algorithms import NumPyMinimumEigensolver\n",
    "\n",
    "numpy_solver = NumPyMinimumEigensolver()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To find the ground state we coul also use the Variational Quantum Eigensolver (VQE) algorithm. The VQE algorithms works by exchanging information between a classical and a quantum computer as depicted in the following figure.\n",
    "\n",
    "<img src=\"aux_files/vqe.png\" width=\"600\">\n",
    "\n",
    "Let's initialize a VQE solver."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from qiskit.providers.aer import StatevectorSimulator\n",
    "from qiskit import Aer\n",
    "from qiskit.utils import QuantumInstance\n",
    "from qiskit_nature.algorithms import VQEUCCFactory\n",
    "\n",
    "quantum_instance = QuantumInstance(backend=Aer.get_backend(\"aer_simulator_statevector\"))\n",
    "vqe_solver = VQEUCCFactory(quantum_instance)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To define the VQE solver one needs two essential elements:\n",
    "\n",
    "1. A variational form: here we use the Unitary Coupled Cluster (UCC) ansatz (see for instance [Physical Review A 98.2 (2018): 022322]). Since it is a chemistry standard, a factory is already available allowing a fast initialization of a VQE with UCC. The default is to use all single and double excitations. However, the excitation type (S, D, SD) as well as other parameters can be selected.\n",
    "2. An initial state: the initial state of the qubits. In the factory used above, the qubits are initialized in the Hartree-Fock (see the electronic structure tutorial) initial state (the qubits corresponding to occupied MOs are $|1\\rangle$ and those corresponding to virtual MOs are $|0\\rangle$.\n",
    "3. The backend: this is the quantum machine on which the right part of the figure above will be performed. Here we ask for the perfect quantum emulator (```aer_simulator_statevector```). \n",
    "\n",
    "One could also use any available ansatz / initial state or even define one's own. For instance,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from qiskit.algorithms import VQE\n",
    "from qiskit.circuit.library import TwoLocal\n",
    "\n",
    "tl_circuit = TwoLocal(\n",
    "    rotation_blocks=[\"h\", \"rx\"],\n",
    "    entanglement_blocks=\"cz\",\n",
    "    entanglement=\"full\",\n",
    "    reps=2,\n",
    "    parameter_prefix=\"y\",\n",
    ")\n",
    "\n",
    "another_solver = VQE(\n",
    "    ansatz=tl_circuit,\n",
    "    quantum_instance=QuantumInstance(Aer.get_backend(\"aer_simulator_statevector\")),\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The calculation and results\n",
    "\n",
    "We are now ready to run the calculation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=== GROUND STATE ENERGY ===\n",
      " \n",
      "* Electronic ground state energy (Hartree): -1.857275030145\n",
      "  - computed part:      -1.857275030145\n",
      "~ Nuclear repulsion energy (Hartree): 0.719968994449\n",
      "> Total ground state energy (Hartree): -1.137306035696\n",
      " \n",
      "=== MEASURED OBSERVABLES ===\n",
      " \n",
      "  0:  # Particles: 2.000 S: 0.000 S^2: 0.000 M: -0.000\n",
      " \n",
      "=== DIPOLE MOMENTS ===\n",
      " \n",
      "~ Nuclear dipole moment (a.u.): [0.0  0.0  1.3889487]\n",
      " \n",
      "  0: \n",
      "  * Electronic dipole moment (a.u.): [0.0  0.0  1.38894899]\n",
      "    - computed part:      [0.0  0.0  1.38894899]\n",
      "  > Dipole moment (a.u.): [0.0  0.0  -0.00000029]  Total: 0.00000029\n",
      "                 (debye): [0.0  0.0  -0.00000075]  Total: 0.00000075\n",
      " \n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.algorithms import GroundStateEigensolver\n",
    "\n",
    "calc = GroundStateEigensolver(qubit_converter, vqe_solver)\n",
    "res = calc.solve(es_problem)\n",
    "\n",
    "print(res)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can compare the VQE results to the NumPy exact solver and see that they match. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=== GROUND STATE ENERGY ===\n",
      " \n",
      "* Electronic ground state energy (Hartree): -1.857275030202\n",
      "  - computed part:      -1.857275030202\n",
      "~ Nuclear repulsion energy (Hartree): 0.719968994449\n",
      "> Total ground state energy (Hartree): -1.137306035753\n",
      " \n",
      "=== MEASURED OBSERVABLES ===\n",
      " \n",
      "  0:  # Particles: 2.000 S: 0.000 S^2: 0.000 M: 0.000\n",
      " \n",
      "=== DIPOLE MOMENTS ===\n",
      " \n",
      "~ Nuclear dipole moment (a.u.): [0.0  0.0  1.3889487]\n",
      " \n",
      "  0: \n",
      "  * Electronic dipole moment (a.u.): [0.0  0.0  1.3889487]\n",
      "    - computed part:      [0.0  0.0  1.3889487]\n",
      "  > Dipole moment (a.u.): [0.0  0.0  0.0]  Total: 0.\n",
      "                 (debye): [0.0  0.0  0.0]  Total: 0.\n",
      " \n"
     ]
    }
   ],
   "source": [
    "calc = GroundStateEigensolver(qubit_converter, numpy_solver)\n",
    "res = calc.solve(es_problem)\n",
    "print(res)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using a filter function"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sometimes the true ground state of the Hamiltonian is not of interest because it lies in a different symmetry sector of the Hilbert space. In this case the NumPy eigensolver can take a filter function to return only the eigenstates with for example the correct number of particles. This is of particular importance in the case of vibrational structure calculations where the true ground state of the Hamiltonian is the vacuum state. A default filter function to check the number of particles is implemented in the different transformations and can be used as"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "=== GROUND STATE ENERGY ===\n",
      " \n",
      "* Vibrational ground state energy (cm^-1): (1e-11+0j)\n",
      "The number of occupied modals is\n",
      "- Mode 0: [0.0, 0.0, 0.0, 0.0]\n",
      "\n",
      "\n",
      "\n",
      "=== GROUND STATE ENERGY ===\n",
      " \n",
      "* Vibrational ground state energy (cm^-1): (2555.376153240112+0j)\n",
      "The number of occupied modals is\n",
      "- Mode 0: [0.9999999999999999, 0.9999999999999999, 0.9999999999999999, 0.9999999999999999]\n"
     ]
    }
   ],
   "source": [
    "from qiskit_nature.drivers.second_quantization import GaussianForcesDriver\n",
    "from qiskit_nature.algorithms import NumPyMinimumEigensolverFactory\n",
    "from qiskit_nature.problems.second_quantization import VibrationalStructureProblem\n",
    "from qiskit_nature.mappers.second_quantization import DirectMapper\n",
    "\n",
    "driver = GaussianForcesDriver(logfile=\"aux_files/CO2_freq_B3LYP_ccpVDZ.log\")\n",
    "\n",
    "vib_problem = VibrationalStructureProblem(driver, num_modals=2, truncation_order=2)\n",
    "\n",
    "qubit_covnerter = QubitConverter(DirectMapper())\n",
    "\n",
    "solver_without_filter = NumPyMinimumEigensolverFactory(use_default_filter_criterion=False)\n",
    "solver_with_filter = NumPyMinimumEigensolverFactory(use_default_filter_criterion=True)\n",
    "\n",
    "gsc_wo = GroundStateEigensolver(qubit_converter, solver_without_filter)\n",
    "result_wo = gsc_wo.solve(vib_problem)\n",
    "\n",
    "gsc_w = GroundStateEigensolver(qubit_converter, solver_with_filter)\n",
    "result_w = gsc_w.solve(vib_problem)\n",
    "\n",
    "print(result_wo)\n",
    "print(\"\\n\\n\")\n",
    "print(result_w)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<h3>Version Information</h3><table><tr><th>Qiskit Software</th><th>Version</th></tr><tr><td><code>qiskit-terra</code></td><td>0.19.0.dev0+105cba3</td></tr><tr><td><code>qiskit-aer</code></td><td>0.9.0</td></tr><tr><td><code>qiskit-ignis</code></td><td>0.7.0.dev0+9201ed8</td></tr><tr><td><code>qiskit-nature</code></td><td>0.2.0</td></tr><tr><td><code>qiskit-finance</code></td><td>0.3.0</td></tr><tr><td><code>qiskit-optimization</code></td><td>0.3.0</td></tr><tr><td><code>qiskit-machine-learning</code></td><td>0.3.0</td></tr><tr><th>System information</th></tr><tr><td>Python</td><td>3.8.10 (default, May 19 2021, 11:01:55) \n",
       "[Clang 10.0.0 ]</td></tr><tr><td>OS</td><td>Darwin</td></tr><tr><td>CPUs</td><td>2</td></tr><tr><td>Memory (Gb)</td><td>12.0</td></tr><tr><td colspan='2'>Mon Jul 26 13:47:51 2021 EDT</td></tr></table>"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<div style='width: 100%; background-color:#d5d9e0;padding-left: 10px; padding-bottom: 10px; padding-right: 10px; padding-top: 5px'><h3>This code is a part of Qiskit</h3><p>&copy; Copyright IBM 2017, 2021.</p><p>This code is licensed under the Apache License, Version 2.0. You may<br>obtain a copy of this license in the LICENSE.txt file in the root directory<br> of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.<p>Any modifications or derivative works of this code must retain this<br>copyright notice, and modified files need to carry a notice indicating<br>that they have been altered from the originals.</p></div>"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import qiskit.tools.jupyter\n",
    "\n",
    "%qiskit_version_table\n",
    "%qiskit_copyright"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
